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SUMMARY  OF  RESEARCH 

Various  aspects  of  the  propagation  of  an  intense  laser  beam 
through  the  atmosphere  are  considered.  The  basic  laser-fluid 
equations  are  presented  and  a  linearized  anal /sis  of  these  equa¬ 
tions  is  given  which  predicts  a  very  low  power  threshold  for 
I  rueckner-Jorna-typc  convective  instabilities.  Another  class 
of  instabilities  is  predicted  to  be  of  more  practical  importance 
than  the  convective  instabilities  and  an  effective  Reynolds  number 
is  derived  which  nay  help  to  characterize  these  turbulent  insta¬ 
bilities.  A  computer  solution  of  the  full  set  of  non-linear 
equations  is  described,  and  the  concept  of  "utility  analysis" 
of  numerical  differencing  schemes  is  introduced.  With  the  com¬ 
putation  scheme  used,  the  laser  pulse  could  be  followed  for  only 
10  seconds;  so  enormous  energy  was  put  into  the  pulse  to  en¬ 
hance  the  interaction  with  the  fluid.  Thus  the  initial  pulse 
distortion  could  be  observed.  Analytical  evaluation  of  the  com¬ 
puter  results  produces  a  detailed  quantitative  check  and  suggests 
tnat  a  combination  of  analytic  and  numerical  methods  would  al¬ 
low  a  pulse  to  be  conveniently  followed  for  much  longer  periods 


Id 


of  time. 
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I  INTRODUCTION 

The  present  paper  and  the  following  paper1  arc  con¬ 
centr'd  with  the  distortions  of  a  laser  beam  produced  by 
density  and  thermal  variations  in  a  fluid  medium.  These 
distortions  have  been  the  subject  >f  many  investigations 
which  can  be  classitied  into  two  groups,  depeuding  ujion 
whether  time  dependence  is  considered.  Most  available 
experiments  are  conveniently  understood  by  refen  nee  to 
theoretical  studies  of  the  gross  effects  of  thermal  dep¬ 
osition  and  fluid  motion  which  assume  that  a  steady 
state  will  be  achieved  for  the  deduction  and  distortion 
of  the  laser  heam.  On  the  other  hand  Brueckner  and 
.forna2  have  discovered  I  hat  some  of  flic  tolutions  for 
beam  projia  gat  ton  arc  unstable  so  that  under  certain 
conditions  a  steady  state  nr>y  not  develop.  The  Brueck- 
ner-.Iorua  instabilities  were  discovered  in  a  linearized 
nnalvsis.  but  the  threshold  for  such  instabilities  has  not 
been  discussed  previously.  The  Brucckiter-.lornn  in¬ 
stabilities  are  not  observed  in  practice  because  they  are 
convective  instabilities  and  cannot  levelop  within  typical 
distances  allowed  for  the  propagation  of  beams  In  the 
laboratory. 

In  the  present  paper  the  linearized  analysis  is  present- 
od.  keeping  the  four-photon  con  ding  Induced  bv  periodic 
fluctuations  in  the  dielec  t  ie  constant,  so  that  among 
other  things  the  thresho!  1  behavior  will  be  exposed.  The 
basic  equations  are  given  in  Sec.  11  and  the  linearized 
analysis  is  presented  in  Sec.  111.  In  Sec.  IV  another 
class  of  instabilities  is  predicted  to  lie  of  more  practical 
Importance  than  the  convective  inst  ibilities  and  an  ef¬ 
fective  Reynolds  number  is  derived  which  may  help  to 
characterize  these  turbulent  Instabilities. 

The  companion  paper1  describes  a  computer  solution  of 
the  full  set  of  nonlinear  laser-fluid  equations.  Knornious 
jiower  was  presumed  for  the  laser  beam  In  order  to 
drive  the  laser-fluid  interaction  as  fast  as  possible  in  an 
attempt  to  watch  the  onset  of  distortions  of  the  beam. 

II  BASIC  EQUATIONS 

When  an  intense  laser  beam  propagates  through  a  fluid, 
manv  inlcrsting  phenomena  take  jilace.  Ibis  laser- fluid 
system  can  be  described  bv  a  macroscopic  model  whi'~i- 
Involves  Maxwell’s  equations,  the  Navier-Stokes  equa - 
tion,  an  energy  conservation  equation,  and  the  continuity 
equation  for  fluid  motion.  These  equations,  which  de¬ 


scribe  Die  behavior  of  intense  electromagnetic  beams 
and  the  associated  sound  and  thermal  fluctuations,  are 
coupled  by  stimulated  Raman  scattering,  electrostric- 
tion,  the  high  frequency  Kerr  effect,  absorption  heat¬ 
ing,  and  the  density  and  temperature  dependence  of  the 
dielectric  constant.  In  this  paper  a  systematic  discus¬ 
sion  is  presented  for  an  intense  laser  beam  propagating 
through  air,  which  lias  a  negligible  Kerr  constant.  If 
the  frequencies  are  outside  the  Raman  scattering  range, 
the  instabilities  are  primarily  caused  by  optical-nc- 
coustic  coupling  of  the  laser  beam  and  the  gases.  These 
effects  are  of  long  dura  lion  compared  lo  those  of  self- 
focusing.  As  the  beam  passes  through  air,  the  Intensity 
profile  induces  a  luiiiuiiifiu  iu  temperature  gradient 
transverse  to  the  prnpnj  itlng  direction  of  flic  beam,  due 
to  the  energy  absorption  from  the  beam.  Tills  thermal 
nonequi librium  and  electrostrictlon  together  cause  the 
generation  of  a  density  gradient  and  hence  a  sound  wave. 
These  density  changes  react  back  oil  the  Inc* dent  beam 
through  changes  In  the  dielectric  constant. 

The  equations  describing  propagation  of  electromagnetic 
radiation  and  the  equations  describing  fluid  behavior  are 
widely  known.’'5  Dropping  unimportant  terms  from  the 
full  equations,  we  take  the  following  set  of  nonlinear 
coupled  partial  differential  equations  for  description  of 
the  macroscopic  representation  of  the  laser-fluid 
system: 

Wave  equation  and  Clnusius-Mosotti  relation: 


vzE=i_i;(fE)  +  i'L(/rE), 

(i) 

f-fL  + 

"W,-  3  ■  wr- 

(2) 

Navier-Stokes  equation: 

Ip! 

/’  —  =  pg-v/M  („  +!ni(, 

(3) 

ki«~»?t'Jv  +  li?  +Tj')V(r-v), 

H) 

15) 

heat  transfer  equation: 
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fluid  continuity  equation: 


(6) 

(7) 


J1  +  v.(pv)  =  0; 

Cl 


(8) 


equation  of  state: 

P  =  P(p,T)  =  RpT.  (9) 

in  the  wave  equation  the  term  involving  rv,  the  linear 
absorption  coefficient,  is  associated  with  a  model  for 
the  absorption  of  electromagnetic  energy  by  the  fiuid." 
The  absorption  coefficient  is  taken  to  be  independent  of 
the  frequency  of  ttie  electromagnetic  radiation,  so  the 
modei  is  not  vaiid  near  the  resonance  lines  of  the 
molecules  in  the  fiuid.  Aiso  the  model  does  not  include 
kinetic  rate  equations,  so  saturation  effects  are  not 
considered.  In  the  present  model  ttie  electric  field  E 
wili  he  damped  by  a  factor  exp(-  \m),  where  z  is  the 
direction  of  propagation  and  the  energy  deposited  in  the 
medium  is  taken  to  be  ntL,  where  lL  is  the 

laser  intensity  in  erg/sec  cm*2. 7 


The  terms  fL  and  f,  are  the  linear  and  nonlinear  permit¬ 
tivity  coefficients,  respectively,  p  is  the  fiuid  mass 
density,  and  v  is  the  velocity  of  a  "material  element"  of 
the  fiuid.  The  convective  derivative  l)  1)1  follows  the 
motion  of  a  material  “particle"  of  the  fiuid  relative  to  a 
fixed  coordinate  system  and  is  expressed  in  the  form 


1)1 


P(J 

?l 


+  v*  v( 


). 


The  vector  g  is  the  gravitational  acceleration  vector  and 
rj  and  rj'  are  the  shear  and  coinpressional  viscosity  co¬ 
efficients,  respectively. 

The  electrostrictive  force  density  f„  Is  given  by 


where  afl  is  the  interaction  stress  tensor  for  the  elec¬ 
tromagnetic  field  and  the  fluid  and  the  angular  brackets 
indicate  a  time  average  over  several  optical  periods.  A 
derivation  of  the  stress  tensor  is  given  on  p.  67  of 
Ref.  4: 


=  J*,,  *'*,*,. 

hut  this  tensor  is  not  strictly  correct  for  optical  fields 
because  an  isothermal  constraint  was  imposed.  A  simi 
lar  derivation  with  an  isentropic  constraint  gives  the 
same  result,  except  that  the  partial  derivative  (?r/i'p)T 
at  constant  temperature  is  replaced  by  the 

derivative  at  constant  entropy.  The  difference  in  these 
two  constraints  is  contained  in  the  thermodynamic 
relation 


(10) 


The  difference  term  in  Eq,  (10)  is  very  small  because, 


for  gases,  (.V/ar),  =r 0  The  term  (yCf/(U'\  -  con¬ 
tributes  a  factor  of  roughly  (.  Actually,  neither  con¬ 
straint  is  strictly  valid.  bnl  corrections  would  he  small 
and  would  necessitate  a  detailed  examination  of  fluid 
boundary  layers  and  Ihc  explicit  mechanisms  of  heat 
deposition  in  the  control  volume. 

The  thermodynamic  quantities  appearing  in  the  above 
equations  are  Cv  and  C ,,  the  specific  heats  in  erg'gdeg 
at  constant  volume  and  pressure,  respectively,  and  y  is 
the  ratio  of  specific  heats,  Ct  Cv\  (i  Is  the  thermal  ex¬ 
pansion  coefficient,  —  (l/p)(?p/dT)p\  r,  Is  the  isentropic 
velocity  of  sound.  [{dP/dp),\l,*\  and  k  is  the  thermal 
conductivity  of  the  fluid 

Due  to  the  complexity  of  the  iaser-fiuid  equations  shown 
above  ,  it  is  not  possible  lo  obtain  exact  solutions  analy¬ 
tically.  The  i  nearized*  solutions  have  been  discussed2 
and  a  number  of  computer  solutions  have  recently  been 
given  by  various  groups."  !n  Sec.  lil  a  linearized  analy¬ 
sis  of  this  set  of  equations  is  presented. 

Ml.  LINEARIZED  ANALYSIS 

Linearized  analysis  is  a  standard  j.  -rturbation  technique. 
In  this  scheme  it  is  assumed  that  each  of  the  dependent 
variables  in  the  problem  can  be  expressed  as  the  sum  of 
its  slowiv  varying  zeroth  order  component  and  a  small 
first-order  correction.9  In  this  way,  a  set  of  linear 
equations  for  small  disturbances  is  obtained.  This  ap¬ 
proach  to  the  analysis  of  the  Iaser-fiuid  system  was 
first  investigated  by  Druerkner  and  Jorna.*  In  the  pres¬ 
ent  approach,  two  variables  are  used  to  describe  the 
perturbed  electromagnetic  field,  one  for  the  component 
of  the  field  which  is  vibrating  in  phase  with  the  primary 
beam  and  one  for  the  component  out  of  phase.  In  this 
way,  the  four -photon  coupling  induced  by  periodic  fluc¬ 
tuations  in  the  dielectric  constant  can  be  included.  This 
coupling  was  not  included  in  the  original  formulation 
given  by  Brueckner  and  Jorna  and  accounts  for  the  ab¬ 
sence  of  a  threshold  ii  their  analysis.  The  Dispersion 
relation  for  these  linearized  equations  has  been  •”"» ‘Mat¬ 
ed  and  is  more  complicated  in  structure  than  that  pre¬ 
sented  by  Brueckner  and  Jorna.  For  propagnticn 
through  air,  however,  the  numerical  differences  are 
minor.  The  wave  with  the  largest  growth  rate,  resulting 
from  resonant  intciactions  between  scattered  electro¬ 
magnetic  waves  and  the  thermal  wave,  propagates  al¬ 
most  perpendicularly  to  the  laser  beam.  The  direction 
is  such  lb:,  t  the  change  in  frequency  of  the  scattered 
electromagnetic  wave  and  the  frequency  of  the  thermal 
wave  (which  is  zero)  are  approximately  the  same. 

Writing 

E  * E(n»  ■*  E<iu 

PBPo  -*  Pi. 

(ID 

rsTo  +  r,, 

f0«?'t(0)  +,2(0l(^(0l)«» 

and  taking 

E,„,  =  ieyE0cxp[i(wj./ - fctz)|exp(- imz)  +c.c. ,  (12) 

E,  =  Jr,{/exp[r(w./  -  k,  •  x)  | 

+  ^exp[-  i(u>J  -k.  »x)|}exp(-  \ oz)  +c.  c. ,  (13) 
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(r)~\(r  )e*p[f'M-k'X)|  U.c,, 


where 

u)L  «  iaser  frequency, 

(15) 

k.*kLr,ik, 

and  treating  the  velocity  v  to  be  a  first-order  quantity, 
we  obtain  iroin  F.qs.  (1)— (9)  the  following  four  equations 
relating  f,  t>,  p\  and  /': 

K'o,  «  i'Hof'M  -  r**!l  f  -*  H'AoSMh: 

■'C2AK0ol)r.'  *(' 0;  (1C) 

2(oi^o a'.)/  +  '  fj(o*A  ~  r2A'2  |.i' 

+  (S AF0ulh>'  >!)T'=  0;  (17) 

(or(f0,),/Jf0|^  +  [,vr(^<),/2F0U'  t  {/(Ct(y  -  l)//?|u>}p' 

+  ( -  kA'2  -  ii>0Cvu')T’  —  0;  (18) 

U('n.  -  IX'o,  4  2)EUA*J|/  4  [i(fo#  -  1  Kr(W  4  2)E0**|tf 
+  {u-2  -  I  .Vwfc'  -  [«2  -  -  1 )  U  O-'}/'* 

•*  {[-«a/3p„  +  ■*('„,  -  l)»E*|fr2}r  -0.  (19) 

Throughout  these  equations,  E„  has  been  written  In 
place  of  IE„I ,  Since  only  the  magnitude  ol  the  complex 
amplitude  appears,  /-.'0  nay  be  considered  to  be  real  and 
positive  without  loss  of  generality,  tor  convenience  in 
writing  these  equations,  the  following  notation  has  been 
introduced: 

«* 1  »•?  V.  (20) 

so  that  ii  is  the  isothermal  speed  of  sound  in  the  un¬ 
perturbed  medium,  and 

(!n 

.Vt(2r,  -i  rj')//»0- 

The  consistency  condition  for  the  linearized  equations 
is  the  vanishing  of  the  determinant  of  the  coefficients  of 
/.  A',  p'.  and  7’'.  Thus  we  obtain  tie  following  re  Li  lion 
between  the  frequency  u>  and  the  wave  vector  k,  the  dis¬ 
persion  relation  for  the  system: 

{(u.  -/«'*%** -/.VW*J  -U2  -A’)A-2|-(y-  D(«2 -«’)cvA-2} 

*(U.  ~  i'llO.^nU.W.2  +  £.U-,2)  '2^t  | 

=  “(A/.,2[Jpu/l(u»  -  in’k2)  ■e/(oc»/0,/3/Ct)(«2  -  /?') J 

-t  /i{[(y-  1)  +i(ocM,1<//>uCv)[ii>2  -f.Vwfr2 

-(«2-zl')/.’2|!)iE2(£.u'!:  s  f.u-.2)^.  (22) 

In  order  to  put  Kq.  (22)  in  the  slightly  more  compact 
form  shown  above,  the  following  additional  abbrevia¬ 
tions  have  been  introduced: 

*'  sK//>riC„, 

hi*  ('o,.)l/J  Index  of  refraction  of  unperturbed 
medium, 


h  power  In  incident  laser  beam/unit  area, 

A'K\('o.-DAFl,  (23) 

W' ■•!('«, 

2wt£.  srJt-2  -f0,w2, 

2  *e2/.'! 

For  coni|Kirison,  tiic  dispersion  relation  obtained  by 
Brueckner  and  dorna2  for  frequencies  outside  the  Raman 
scattering  range  is 

u>(u>2-r2fr2-/.Vu>fcJ)£t  -t/t(i/w +/«)**=  0,  (24) 

where 

Cz  b  l**  “  (2^to,/«o.r)/tfr2t  |  -  4(Mn,fct/e*)»[w  -  <r/«J*t|2, 
v*p»A2k\  <c,  (25) 

6s2A/t/'2oA-2/MJ,C„. 

Although  Fq.  (22)  is  considerably  more  complicated  in 
structure  than  Kq.  (24),  the  general  features  of  the  two 
equations  are  the  same  As  a  lirst  approach  to  the 
analysis  of  (22),  one  should  realize  that  the  power  in 
the  primary  laser  beam  is  proportional  to  Ej.  Thus, 
tiie  free  modes  of  the  system  can  lie  obtained  by  ietting 
i-o~* 0.  With  no  power  in  the  incident  beam,  therefore, 
(22)  reduces  to 

{[w  -  iU/p0C,)fca  |(u>2  -  /ATwfc2  -  nV)  -  (y  -  1  )//2u»A’2} 

x[r2(A’2  4A2  4  2kLk,)-fL{ot(u>L  4u>)2| 

x[cJ{k2L  4A-2-2AtA',)-fUo)(a,i-w)2|  =  0.  (2G) 

The  first  factor  in  (2(1)  contains  a  nonpropagating  ther¬ 
mal  wave  and  two  damped  sound  waves  coupled  by  the 
term  (y  -  1  )n2JL>k2.  The  last  two  factors  correspond  to 
the  four  free  inodes  for  scattered  undamped  electro¬ 
magnetic  waves: 

w/wt  ffi[l  4*2  A’i  4-2ff<f,/|*tll,/2,  (27) 

where  o=  1 1,  specified  by  given  values  of  k1  and  A-,. 

Two  of  the  roots  are  low  frequency  (ui«wt),  whereas 
the  other  two  have  frequencies  of  the  same  order  as  the 
laser  frequency,  it  is  clear  Unit  the  roots  at  the  high 
frequency  should  he  eliminated,  because  it  has  boon 
assumed  previously,  in  evaluating  time  averages,  that 
the  perturbed  solutions  vary  much  more  siowiy  than  the 
optical  waves.  Therefore,  factors  i ike 

^■’-'moiW.2  (28) 

will  be  replaced  by 

^uoi^u^Hr^A’,  -  a’,),  (29) 

where 

<"l  * <'/[*' tt oi  I* / 2  velocity  of  light  in  the  medium. 

Thus  the  free  inodes  wili  include  three  thermal -sound 
waves  and  two  electromagnetic  waves. 

All  of  the  terms  in  (20)  resu't  from  the  ieft-hand  side  of 
(22)  because  the  right-hand  side  is  proportional  to  the 
power  in  the  primary  laser  beam.  Now,  as  the  power  in 
Uie  laser  beam  is  turned  on,  the  right-hand  side  couples 
the  live  free  modes  described  by  (20).  Additional  tiny 
coupling  arises  inside  the  left-hand  side  itsilf  through 
the  A',  If,  and  r3(0l  terms. 
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For  detailed  consideration,  the  case  of  a  primary  laser 
beam  at  10.  C  p  propagating  through  air  at  approximate¬ 
ly  10  C  and  at  standard  pressure  will  be  discussed.  The 
numerical  values  for  the  parameters  appearing  In  the 
dispersion  relation  are10 

uit=I  773x10“  sec*1, 

kL  —5920  cm*', 

^=1.25  xio-1  g  'em1, 

,V=  0.2U4  cm1  sec, 

/?  =  7. 07  x  10*1  deg*1, 

C„=  7.  M3  xlO"  erg/gdog, 

k'  s=0.  28  cm2/scc, 

n0=  I  +2.82x10“, 

u2  -8. 39  xlO"  (cm  sec)2, 

y-  1.4, 

rt(0)=l  45.05  xl0*\  (3°: 

rr  =  (3  x  10*7  cm*')«0, 
nrjOrcu20/Cv=  3.  831  xlO’ciy,, 

"70,/ftA  =0.988ryo, 

‘rJ(0,/-2  =  (I.37xIO*2!  sccVgcm2)/,. 

A  -  0.452  cmVg  ■*  (1-  10x10*”  cm  sec1  g2)/,, 

A(y  -  1 )  2/35=24. 7  cm1  deg  g 

^(G.OOxlo*21  cm  sec1  deg  g1)/ L, 

1  4  5. 05  xio*1  -♦  (1 , 37  x  10*21  sec’/g  cm2)/,, 

A •  (2. 84  x  10*ls  cm2  sec  g)/t  4  (l .  38  x  10*ln  secVg1)//, 

B  ~'nn>£l/To  =  “(•*• 84  scc’/g  cni2dcg)/t, 

B'  =5  -(0.  602  xlO*1"  secYgtfi. 


instead  of  A-,,  and  defining 

T*Au>Lflntl,JJfrCp,  (33) 

which  corresponds  to  the  power  parameter  used  in  Hef. 

2,  the  dispersion  relation  can  be  written  in  the  form 

(uj  -  1  k'A'2)(u>2  -  /'  Vo)/;2  -  it7k2)  -  (y  -  1  )i'?u.'A'2 

=■  J»Ti'JAr*{[tat-  \incL-rL{hv  -t  It2/ 2A*t)|*1 

-[u>  -  \inrL-rL(hv  -  k3/2kL )|*'J.  (34) 

The  problem  at  this  stage  is  the  determination  of  the 
maximum  growth  rate  of  any  Fourier  component  of  a 
distortion  of  Die  plane  wave  as  a  function  of  the  ab¬ 
sorbed  power  from  the  beam.  That  is,  one  must 
solve  the  dispersion  relation  for  the  frequencies  ns  a 
function  of  A>,  c,  lL,  and  the  characteristic  parameters 
of  the  medium,  ard  then  find  the  maximum  value  of 
-lino  for  real  v  and  k  with  IH  <1.  Such  a  problem  can¬ 
not  be  solved  analytically  without  further  approxima¬ 
tions.  One  region  of  interest  would  be  the  high  power 
limit,  where  the  driving  term  would  overwhelm  the 
losses  resulting  from  thermal  conduction,  viscosity, 
and  the  absorption  of  electromagnetic  energy.  In  that 
case,  all  the  imaginary  terms  in  Kq.  (34),  with  the  sole 
exception  of  the  i  immediately  preceding  r,  can  be 
dropped,  reducing  the  dispersion  relation  to  the  (orm 
employed  by  Bruecker  and  Jorna  in  Kq.  (45)  of  Hef.  2. 

To  proceed  analytically,  Brueekner  and  Jorna  neglected 
the  term  in  the  second  set  of  brackets  on  the  right-hand 
side  of  (34),  and  assumed  that  the  maximum  growth  rate 
would  occur  somewher  on  the  curve  in  the  v,  A-  plane 
determined  bv  the  constraint 

He[u.>  -  cL(hv  +  A,J  '2/.'t)|r-0.  (35) 

Along  that  curve,  the  maximum  growth  rate  is 

(-Imw)ml,=  J/f(1.08),  (36) 


In  the  above  list  a  dimensionless  absorption  constant  n0 
of  order  unity  has  been  introduced  and  the  power  /,  is 
in  units  of  erg  sec  per  cm2.  Now,  using  these  numerical 
values,  one  finds  Dint  the  power -dependent  terms  are 
very  small  for  power  fluxes  less  than  10  MW/cm2,  ex¬ 
cept  for  Uie  term  which  represents  energy  absorption. 

In  other  places  in  the  dispersion  relation,  the  power- 
dependent  terms  are  connected  with  the  nonlinear  index, 
and  will  be  omitted  in  the  following.  (The  terms  omitted 
are  related  to  .:clf -focusing  in  .1  manner  described  by 
Brueekner  and  Jorna.2)  This  neglect  of  the  nonlinear 
index  and  of  the  weak  dependence  of  the  optical  coeffi¬ 
cients  of  gases  on  the  temperature  for  fixed  density 
allows  tlie  simplification  of  Uie  dispersion  relation  to 

[(u-  -  iK’lSlU  -  i.Vu- A-2  -  I'V)  -  (y  -  1  )if2wA’2 1 

/\cL A’.  -  -  u>.) 

=  -{Ab2uLIL  2fv,,)[(rtA-.  -  w.)  +{cLk.  -u*.)| 

Xj^po'Hw  -iK’b2)  4  jnrn0(hi2/Cv\.  (3 i ) 

In  addition,  for  air  at  reasonable  powers,  Uie  term 
J/i,  /HiL'-irT2)  on  the  right-hand  side  is  negligible  com¬ 
pared  to  'yen,, (hi2  O',,.  Introducing  the  variable 

v  *  A-j/Ar,  (32) 


which  corresponds  to 

(.T),/,/°.g751_S\/  1  \l/2 

'  \0. 56305/\0. 93003 /  ' 


(37) 


(These  results  differ  from  those  in  Hef.  2,  which  are 
erroneous.)  There  is  no  assurance  that  the  ariual  maxi¬ 
mum  growth  rate  does  lie  along  Die  one-dimensional 
subset  of  the  v,h  plane  assumed  in  Ref.  2.  We  have  con¬ 
ducted  a  search  along  the  line 


1*  4  It  '2hL  0. 

However,  the  result  for  the  maximum, 

|/r(1.00),  (38) 

is  2',  smaller.  No  other  curve-  in  Hie  v,k  plain  has  been 
found  which  allows  an  analytical  search.  Nevertheless, 
one  suspects  that  these  answers  are  sufficiently  close 
and  Unit  further  analytical  effort  Is  not  justltied,  be¬ 
cause  of  tin*  previous  approximations. 

An  interesting  unknown  not  discussed  previously  is  Uie 
power  flux  required  to  stimulate  these  instabililies. 

This  threshold  power  is  c  learly  a  critical  function  of  the 
losses  in  the  system,  which  therefore  renders  it  im¬ 
portant  to  treat  them  carefully.  If  the  second  term  on 
the  right-hand  side  of  Kq.  (34)  is  dropped,  the  instnbil- 


ity  appears  to  have  no  threshold,  because  the  comiuction 
loss,  which  must  be  overcome,  vanishes  as  k  —  0.  iiow- 
e\.  r,  as  *■  —  0,  llie  Stokes  and  anti-Stokes  terms  on  the 
right-hand  side  of  Ftp  (3-1)  tend  to  cancel  each  other, 
and,  therefore,  there  is  a  threshold  (lower  flux  lor 
these  stimulated  thermal  Rayleigh  scattering  instabil- 
lties.  Using  (3-1),  a  computer  search  for  this  threshold 
was  performed  and  led  to  the  result. 

(/e)thr.,f,oid-';0-32y  m W/ cm2.  (3‘J) 

This  threshold  was  located  at  /,■  0,04  cm'1  with  v 

j  1  1  «'10"  .  The  degeneracy  in  the  value  of  v  occurs 
because  of  the  symmetry  property  contained  in  (34), 
r  "  -  t'  implies  iw  —  (lad*.  The  dispersion  relation  of 
iiru<  ckner  and  Jorna  contained  no  threshold  for  the 
convective  instabilities  and,  theroluro,  suggested  that 
inslaiiililies  might  he  present  for  extremely  low  beam 
intensities.  The  J-niW  per  cm2  threshoid  obtained  in  the 
present  analysis  is  cert.iinly  small  in  relation  to  i  it  - 
tensities  available  for  experiments. 

Tin*  presence  of  a  wind  does  not  alter  ti  e  growth  rates 
for  distortions.  This  can  be  easily  seen  by  considering 
the  problem  from  a  frame  of  reference  moving  with  the 
fluid.  A  uniform  beam  remains  a  uniform  beam  in  the 
moving  frame,  although  its  direction  of  propagation  is 
.shifted.  Tills  shill  in  direction  iias  no  effect  upon  the 
stability  discussion.  The  very  iow  convective  instability 
threshold  does  not  mean  that  such  effects  are  easy  to 
observe  or  ..re  of  practical  importance.  An  extremely 
long,  and  carefully  protected,  optical  path  length  would 
be  required  in  order  to  see  these  instabilities  grow  to 
observable  sizes.  Indeed,  such  instabilities  have  never 
been  observed  in  the  laser-fluid  system. 

There  is,  however,  another  class  of  instabilities  in¬ 
volving  the  ha  lance  between  Inertial  terms  and  bouyancy 
forces  in  the  fluid  Such  instabilities  are  dis'-ussed  In 
Sec.  IV. 

IV.  BEAM  FED  TURBULENCE 

Many  experiments  have  been  reported  for  which  the¬ 
ories  assuming  steady-state  beam  profiles,  after  initial 
transients  die  out,  provide  rather  good  explanations  of 
the  principal  features.  However,  ’.hat  is  probably  true 
only  because  lliose  experiments  arc  conducted  at  rcla- 
tiveiy  low  power  fluxes.  Theoretically,  one  expects  a 
time -dependent  state  of  the  system  because  of  the  In¬ 
stabilities  discussed  earlier.  Such  instabilities  are  not 
observed  in  practice  because  they  cannot  develop  within 
tin*  distances  allowed  for  Ui  :  propagation  of  beams. 
However,  alternate  considerations  lor  a  beam  of  finite 
cross  section  suggest  that  the  beam  may  drive  the  fluid 
into  a  time-dependent,  or  turbulent,  state  at  powers 
which  art  noi  completely  unreasonable. 

It  may  be  impossible  to  prove  analytically  Unit  such  a 
turbulent  state  develops,  because  the  investigation  of 
hydrodynamic  stability  is  very  difficult  even  for  Die 
simplest  flows.  However,  an  argument  can  be  made 
from  dimensional  considerations,  an  approach  that 
promises  to  lie  very  useful.  Naim  ly,  lor  the  problem 
of  a  lx  am  of  radius  <i  ami  power  lliix  I  passing,  through 
an,  it  is  pissible  In  estnu.iti  a  parameter  It’,  whieli 
plays  tlic  role  of  an  efiective  Reynolds  number  for  uur 


problem.  It  will  be  shown  that  the  parameter  IF  Lakes  on 
values  of  the  order  of  30000  for  a  beam  with  Intensity 
/•  1  k\V  cm2  of  radius  1  m.  Since  it  Is  known  that  some 
flows  with  Reynolds  numbers  .subsLintlally  lower  Hum 
30000  are  turbulent,  the  flows  for  the  laser-heated 
atmospheric  jxith  should  also  be  expected  to  show  signi¬ 
ficant  time  dependence. 

Consider  the  equations  of  motion  for  the  air  and  the 
equation  governing  heat  transfer,  which  take  the  follow¬ 
ing  form  if  it  is  assumed  that  the  air  can  be  assumed 
incompressible  (that  amounts  to  dropping  terms  of  order 
x1  c2,  where  i<  Is  a  typical  tiow  speed,  and  t ,  is  the 
speed  of  sound;  for  tl  e  problems  under  consideration, 
a  i  2  will  be  less  than  1 0'1,  and  the  Incompressible 
fluid  approximation  will  be  quite  good): 

p(v*  v)v  -  -  V/>-  /37‘,pg  +(p  /p0)i|V2v 
+  (f*/PiiHrj  v), 

p(v.v)=  -(v.  V),,  =-/3p(v.  v)r„  (40) 

pt',(v.  v) r,  -  <y/4  r.U-rr,) 

+  (p  +T)'(V.  v)2|. 

Assuming  the  Reynolds  number  is  high,  the  inertial 
terms  will  dominate  the  viscous  terms  in  the  Nnvier- 
Stokcs  equation.  Thus,  there  must  be  a  balance  between 
the  inertial  terms  and  the  bouyancy  forces,  which  im¬ 
plies  that  (uc  (i  ~n^l  where  p  is  the  density  of  air, 

/iis  the  coefficient  ot  thermal  expansion,  T,  is  a  typical 
value  lor  the  temperature  rise,  and  g  is  the  accelera¬ 
tion  due  to  gravity.  The  pressure  variation  will  be  of 
order  pi/2.  In  the  heat-transfer  equation,  the  convection 
term  will  dominate  the  conduction  term,  and  the  beam 
heating  will  overwhelm  the  viscous  dissipation,  so  that 
•here  must  be  a  balance  between  heat  deposition  from 
the  beam  and  convective  heat  transfer,  which  implies 
that  pC^nT,  n~nl.  Combining  these  two  relations  we 
find  that  i/1  -  nfri'lg/pCf,.  With  tills  expression  for  w,  we 
then  define  a  parameter  if,  which  is  expected  to  in¬ 
dicate  regimes  where  steady  flow  and  where  tinic-dc- 
pendent  flow  may  be  anticipated.  IF  is  an  estimate  uf  the 
relative  importance  of  iuerti.il  terms  to  viscous  terms 
in  controlling  the  How: 

W  =  a  it p0/r) 

—  n(po/T)Ko/3f72/^/pC^)* /5.  (41) 

Fora  beam  with  /=  1010  erR./cni*  sec,  n  =  100  cm,  T0 
-273  K,  and  a  -  3  *1  O'"  cm'1,  11’^  30000. 

For  many  experiments  described  in  the  literature,  the 
values  of  IF  are  much  smaller,  and  thus:  one  would  not 
expect  any  turbulent  iluid  flow  to  be  observed.  For  ex¬ 
ample,  in  the  original  experiment  of  Gordon  Leite, 
Moore,  Porto,  and  Whinnery11  the  parameter  It’  Likes 
on  a  value  about  10'1,  and  in  the  more  recent  experi¬ 
ments  of  Smith  and  Gebhardt,12  IF  is  of  order  10. 

The  parameter  IF  introduced  here  is  different  from  the 
Grasshof  number,  which  is  referred  to  in  some  discus¬ 
sions  of  the  convective  flows  set  up  by  the  absorption  of 
energy  from  a  laser  beam. 13  In  fact,  the  court  piu.il 
basis  lor  using  the  Grasshof  number  in  a  discussion  at¬ 
tempting  to  explain  the  transition  between  smooth  How 
and  time -dependent  flow  seems  less  relevant  because 
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tho  Grnsshof  number  appears  to  be  more  sensible  when 
trie  thermal  bouynnt  forces  are  balanced  try  viscous 
forces.  In  the  present  discussions,  the  thermal  bouynnt 
forces  are  balanced  bv  inertial  effects.  It  turns  out  that 
the  number  W’  is  essentially  the  square  root  of  the 
Grasshof  number. 14 

We  are  planning  experiments  to  determine  the  critical 
value  of  U’,  ll'rr,  which  determines  the  onset  of  turbulent 
convective  fiows  for  the  geometry  appropriate  to  laser 
beam  transmission.  It  is  also  our  aim  to  attempt  a  the¬ 
oretical  evaluation  of  this  critical  vaiue.  At  the  present 
time  we  can  oniy  speculate  that  H'cr  may  be  between  101 
and  10V*  The  theoretical  approach  appears  fairly  dif¬ 
ficult  because  the  question  of  the  stability  of  flows  even 
without  heat  sources  has  only  been  answered  theoreti¬ 
cally  for  very  simple  geometries. ,f,,IT  The  question  of 
stability  for  fluids  which  are  heated  or  cooled  appears 
to  have  been  treated  mainiy  for  cases  in  which  the  fluid 
would  be  motionless,  and  has  not  been  explored  for  a 
problem  like  the  present  one. n’1"  The  first  part  of  that 
problem  would  be  to  determine  a  steady-state  fiow  pat¬ 
tern  for  a  fluid  with  a  distributed  heat  source  within  a 
right  circular  cylinder  with  its  axis  aligned  at  some 
angle  to  the  vertical.  For  the  case  of  a  horizontal  cyl¬ 
inder  of  Infinitely  smali  radius,  the  flow  pattern  has 
been  calculated  by  Yih.le,,° 

Unfortunately,  however,  that  solution  is  not  of  great 
vaiue  for  the  present  problem  because  the  size  of  the 
cylinder  radius  is  a  critical  parameter.  Nevertheless, 
it  is  expected  that  Yih’s  solution  will  assist  in  obtaining 
the  asymptotic  form  of  the  steady-state  fiow  at  large 
distances  from  the  laser  beam  cylinder.  Once  that 
time -independent  flow  pattern  has  been  determined,  the 
linearization  of  the  hydrodynamic  equations  for  pertur¬ 
bations  from  the  flow  pattern  will  lead  to  an  eigenvalue 
problem,  which  eventually  will  yieid  a  critical  vaiue  for 
W.  Ostracii21  suggests  that  the  eigenvalue  problem  can 
be  bypassed  as  the  stability  of  fullv  developed  natural 
convection  flows  can  be  found  by  using  the  appropriate 
velocity  profile  in  the  classical  theory  of  hydrodynamic 
stability.  This  assertion  rests  upon  his  analysis  of  the 
stability  of  free  convection  above  a  fiat  iieated  plate, 
where  instability  first  appears  for  a  Reynolds  number 
of  283. 

Above  the  threshold  for  beam-induced  turbulence,  gov¬ 
erned  by  Wcr,  general  arguments"  lead  to  a  size  for  the 
smallest  eddies,  a(Wct/W)ili.  Fora  1-m  beam,  if  lFcr 
should  be  about  10’,  then  the  eddies  might  have  sizes 
as  small  as  7  cm  for  a  power  flux  of  1  kW/cm1,  The  as¬ 
sociated  density  fluctuation  would  then  be  expected  to 
result  in  considerably  increased  scattering  of  the  beam. 

The  arguments  presented  here  show  that  there  are  sub¬ 
stantially  more  important  sources  of  instability  in  tiie 
laser-fluid  system  than  those  discussed  In  earlier 
linearized  analyses.  It  is  felt  that  these  fluid  instabil¬ 
ities  will  be  enhanced  bv  their  interaction  with  the  scat¬ 
tering  of  the  laser  oeam,  because  of  the  general  result 
that  instabilities  in  fluids  result  if  the  heating  of  the 


fluid  is  greater  in  those  regions  where  the  density  of  the 
fluid  is  greater. 21 

At  the  present  time  we  can  only  outline  the  general  na¬ 
ture  of  the  effects  to  be  expected  above  a  critical  power 
level,  Much  additional  work  clearly  needs  to  be  done, 
both  of  an  experimental  and  theoretical  nature. 
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I  INTRODUCTION 

TI10  work  described  in  this  paper  was  motivated  by  the 
discovery  of  instabilities  in  the  system  of  equations  de¬ 
scribing  electromagnetic  wave  propagation  and  fluid  dy¬ 
namics.  These  instabilities  are  described  in  the  preced¬ 
ing  paper.  A  numerical  solution  of  the  full  set  of  mac¬ 
roscopic  laser -fluid  equations  (without  saturation  ef¬ 
fects)  Is  presented  for  a  Gaussian  laser  pulse  of  enor¬ 
mous  energy.  The  computation  scheme  employed  was 
not  optimal  and  the  pulse  could  be  followed  for  oniv 
lO-’  sec.  Consequently,  enormous  energy  was  pit  into 
the  pulse  to  enhance  the  laser-fluid  interaction  and  drive 
the  onset  of  strong  distortions  and  instabilities. 

Other  computer  solutions  of  the  laser-fluid  equations 
have  recently  been  given'  bv  various  groups  for  physi¬ 
cally  reasonable  psiwers.  Such  calculations  have  typical¬ 
ly  dealt  with  steady -state  beams  without  saturation  ef¬ 
fects.  bu*  h,,ve  included  the  effects  of  gravity.  More  re¬ 
cently  s. dural  ion  effects  have  been  added  and  pulse 
shape  is  being  considered. 

I11  the  course  of  Ibese  studies  it  become  apparent  that 
lliere  was  some  merit  to  introducing  a  new  concept  to 
judge  the  value  of  an  algorithm  for  computing  the  solu¬ 
tions  of  a  system  of  partial  differential  equations.  This 
concept  was  called  "utililv".  and  will  lie  discussed  only 
briefly  in  Sec.  fl.  The  full  details  will  lie  published 
elsewhere.1  The  advantage  of  this  concept  is  that  it  Is 
relatively  easy  to  apply  to  complicated  systems  of  par¬ 
tial  differential  equations,  whereas  the  stability  concept 
leads  to  a  very  complicated  procedure  for  deciding  on 
the  value  of  a  numerical  roudne. 

Speed  and  memory  size  in  a  computer  place  certain  re¬ 
strictions  on  one’s  ability  to  investigate  phenomena  in 
the  laser-beam  problem.  In  tin;  attempt  to  calculate  dis¬ 
tortions  of  tlm  type  predicted  by  the  linearized  lnstabili- 
ly  analysis,  cylindrical  symmetry  was  imposed  on  the 
problem  in  order  to  facilitate  the  computer  calculation. 
Had  this  not  been  necessary,  or  had  some  other  inde¬ 
pendent  variable  been  eliminated  rather  than  the  angle 
about  tire  beam  axis,  much  more  pronounced  evidence  of 
beam  and  ftuid  instabilities  would  likely  have  been  ob¬ 
served  for  substantially  lower  powers,  jrowers  that  may 
b<*  achievable.  Arguments  supporting  this  proposition 
arc  contained  in  the  preceding  paper, ' 


In  Secs.  Ill  and  V  we  present  the  results  of  a  calculation 
of  beam  distortion  for  a  very  high  Intensity  pulse  prop¬ 
agating  through  air  for  several  kllomters.  Analytical 
agrumcnls  are  advanced  In  Sec.  IV  which  suggest  that 
the  qualitative  features  of  the  distortions  are  correct, 
which  lends  credence  to  the  computer  output.  Then, 
using  the  computer  results  as  a  check,  the  analytical 
procedure  Is  shown  to  be  adequate  for  detailed  quantita¬ 
tive  calculations.  Rv  combining  the  analytical  procedure 
developed  with  approy* date  computer  support,  the  pulse 
couid  be  followed  much  longer  than  the  10''  sec  de¬ 
scribed  here. 

II.  DESCRIPTION  OF  THE  PROCEDURE 

The  laser-fluid  equations  were  solved  in  the  near-fleld 
region  of  a  laser  pulse,  initially  Gaussian  In  both  r  and 
z,  propagating  through  air  at  1  atm  of  pressure  and  at 
10'C.  A  cylindrical  geometry  was  used  and  cylindrical 
symmetry  ,no  dependence  on  the  ant  le  <?>)  was  preserved 
at  the  price  of  dropping  the  gravity  '  <rm  in  the  Navicr- 
Stukes  equation.  Waving  evlindrir.il  symmetry  amounts 
to  a  considerable  simplification  in  fie  problem,  so  that 
the  Inclusion  of  Hie  free  convection  iffccts  due  to  gravi¬ 
ty  was  nut  attempted  in  this  analysis.  The  problem  de¬ 
scribed  above  amounts  to  a  mixed  initial-boundary  prob¬ 
lem  The  initial  configurations  of  the  laser  beam  and  the 
fluid  are  specified  subject  to  certain  boundary  conditions 
at  !•=■()  which  must  be  satisfied  as  all  times.  Further¬ 
more,  trie  boundary  condition  at  *  =  0  is  time  dependent, 
because  Hie  tall  of  the  Gaussian  must  be  fed  into  the 
spatial  region.  Per  the  numerical  solution,  a  spatial 
mesh  of  grid  points  or  stations  is  useo  to  represent  the 
v?  plane.  At  a  given  instant  in  time,  the  values  of  the 
various  dependent  variables  are  obtained  at  all  of  the 
stations.  The  diffei  once  equations  are  then  employed 
with  these  values  of  the  dependent  variables  to  advance 
a  step  in  time.  This  procedure  is  repeated  over  and 
over  until  the  desired  time  interval  has  been  traversed 
An  explicit  difference  scheme  was  used  In  this  calcula¬ 
tion,  because  such  schemes  arc  simplest  to  handle. 

The  two  major  difficulties  in  using  numerical  techniques 
to  solve  differential  equations  bv  computer  are  error 
growth  and  excessive  compulation  time.  In  order  to  con¬ 
trol  the  error  growth,  the  utility  criterion  mentioned  in 
Sec.  I  has  been  used.  Furthermore,  highly  accurate 
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I'll!.  I.  On-axls  inar’iltufk1  of  the  cifctrlc 
fluid  ifl  is  shown  as  a  function  of  |*-if| 
at  four  different  times.  As  shown  in  Vq. 

(1C) ,  ze  »  rt  loc  ates  the  "center”  of  the 
pulse.  The  solid  curves  show  the  leading 
edge  of  the  pulse  and  the  dotted  curves  de¬ 
pict  the  trailing  edge.  Distances  alone  the 
t  axis  arc  expressed  In  units  of  the  grid 
size:  a*  0.45  km.  The  quantity  l£|  she  wn 
is  defined  in  Kq.  (13) ,  so  that  the  exponen¬ 
tial  damping  factor  is  not  Included  In  ttic 
graphs. 


U-»el 


A2 


1  <AZ *0  4  5  km  )  • 


seven -point  difference  quotient  representations  of  the 
differer.tici  operators  wore  employed  to  reduce  trunca¬ 
tion  error.  In  order  to  handle  the  economic  problem  of 
large  computation  time,  a  certain  amount  of  efficiency 
is  introduced  by  minimizing  the  amount  of  core  storage 
required  of  the  computer.  This  was  accomplished  in 
part  by  using  overiaving  techniques  to  store  several 
pieces  of  information  at  the  same  site  in  ihc  computer. 
Thus  information  is  stored  only  as  long  as  it  is  needed 
and  then  is  replaced  with  current  material  The  com¬ 
putation  time  was  also  minimized  bv  making  us<  of  a 
nonuniform  grid.  Tile  seven-point  difference  relations 
allowed  a  relatively  large  grid  size  without  undue  trun¬ 
cation  error  and  the  nonuni  form  grid  spicing  permitted  a 
greater  grid  density  in  the  region  of  special  interest. 
Thus  an  accurate  solution  could  be  obtained  with  a  min¬ 
imum  of  computation. 

The  laser-fluid  equations  are  given  in  l-.qs.  (1)— (9)  of 
Ref.  1.  As  mentioned  above,  the  gravity  term  wz.4 
droppeu.  Also,  the  thermal  conductivity  «  was  taken  to 
be  constant  because  its  derivatives  are  very  small.  The 
equation  of  state  was  taken  io  he  the  ideal  gas  law  The 
numerical  values  used  for  the  various  parameters  are 
the  sama  as  those  given  for  the  linearized  analysis  in 
Ref.  1,  because  the  same  temperature  and  pressure 
were  used  for  the  undisturbed  medium.  The  laser  fre¬ 
quency  u)L  and  the  dimensionless  absorption  constant  o„ 
were  chosen  to  be  1.773*10M  sec-'  ami  10,  respec¬ 
tively.  The  wave  equation  for  ilnenrly  polarized  light  in 
an  absorptive  medium  Is  reduced  to  the  scalar  equation 

!•:).  (1) 

This  equation  is  an  approximate  equation  describing  an 
electric  field  which  is  polarized  linearly.  Strictly 
speaking,  of  course,  Maxwell’s  equations  do  not  allow 
cylindrlcally  symmetric  linearly  polarized  beams  in 
charge-free  space. 

The  solution  of  (1)  Is  taken  to  be  in  the  form 

£=  i(E’,  +i£j)expi(ufc(-.*t,r)exp(-!nr*)  +  c.c.  ,  (2) 

where  E ,  and  E,  are  slowly  varying  'ancuons  of  r  a..J  , 
and  the  laser  frequency  and  wave  number  are  related  by 


-c*kL  ,  (3) 

where  <Bf  is  defined  in  Kq.  (11 )  of  Ref.  1.  Equation  (2) 
was  substituted  into  (1)  and  the  second  derivatives  of  /:, 
and  Ej  with  respect  to  time  were  dropped.  From  the 
real  and  imaginary  parts  of  til*-  resulting  equation,  cou¬ 
pled  simultaneous  differential  equations  for  E,  and  / 
were  obtained  These  equations  are  also  coupled  to  the 
differential  equations  for  /»,  T,  and  lhc  cylindrical  com- 
jxments  vr  and  r  r  of  the  fluid  velnciiy.  The  use  of  cyiii  - 
drical  symmetry  allows  the  slmplllication  v„^0. 

The  easiest  way  to  obtain  a  variable  grid  size  is  to  in¬ 
troduce  a  transformation  to  a  now  indpemient  variable, 
Thus,  in  order  to  have  more  grid  |»ints  in  the  region  of 
special  Interest,  small  r,  me  nonlinear  transformation 

r/r„-=x/(\ -x)  (4) 

was  employed.  The  scale  vniue  was  chosen  according 
to  the  dictates  of  convenience  and  will  be  shown  in  Ka 
0). 

Because  of  the  symmetry  of  the  problem  and  the  regu¬ 
larity  of  the  differential  equations,  the  following  bound¬ 
ary  conditions  must  be  satisfied  at  v~0: 

?liy  Op  Or,  ?vr  _P:rr  _ . 

Ov  «"» v  ?v  0*  Or  Bz'iz1’1'  0 


The  laser-fluid  differeniial  equations  must,  of  course, 
be  converted  to  difference  equations  before  a  computer 
solution  can  be  attempted.  As  indicated  earlier,  seven- 
point  difference  quotients  were  used  to  represent  the 
differential  operators,  but  these  difference  quotients 
will  not  be  presented.  (The  scheme  was  a  straightfor¬ 
ward  explicit  difference  representation. ) 

Once  the  difference  equations  have  been  written,  It  Is 
vital  to  have  some  criterion  to  determine  useful  time 
step  sizes  and  corresponding  grid  spneings.  It  is  not, 
however,  necessary  to  demand  that  the  difference 
scheme  be  stable  in  the  classical  sense.  In  fact,  the 
classical  notion  of  stability  is  irrcvelant  to  computer 
solut  a  of  inferential  equations,  particularly  for  non- 
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l/rA/>  l2/kL(raMf  +  B/te.  (8) 

Since  Z  Is  given  oniy  approximately  by  (7),  the  laser 
l>ower  4  floes  not  appear.  The  differencing  scheme  em¬ 
ployed  was  not  optimal  and  its  utility  region  Is  entirely 
dominated  by  tiie  eleclric  fleid  development  and  Is  In¬ 
sensitive  to  the  fluid  parameters.  This  was,  however, 
the  scheme  employed  In  the  present  solution. 

Condition  (8)  Is  such  a  strong  constraint  that  one  imme¬ 
diately  wonders  if  it  is  really  necessary  to  obey  it.  Part 
of  the  utility  philosophy  Is  to  obtain  a  constraint  iike  (8) 
and  try  it  on  the  computer.  Then  one  can  try  to  violate 
the  condition,  using  a  larger  A/,  and  compare  the  re¬ 
sults.  This  was  done  for  the  problem  discussed  in  this 
paper  and  no  escape  from  (8)  was  possible  for  the  dif¬ 
ference  scheme  used.  In  fact,  if  the  criterion  was  vio¬ 
lated  by  a  factor  of  order  5  In  A/,  then  classic  Instabil¬ 
ity  phenomena  were  observed  In  the  computer  output. 
Thus,  by  a  stroke  of  bad  iuck,  it  appears  that  (8)  must 
be  obeyed. 

In  order  to  emphasize  the  Implications  of  (8)  for  the 
study  of  the  propagation  of  laser  puiscs,  a  description 
of  the  accessible  parameter  regime  wiii  now  be  given. 


KKJ,  Fur  detailed  comparison,  the  z  profile  of  the  on-axis 
electric  Held  I  1  is  shown  :il  several  times,  'I he  unit  A*  is 
used  lor  distances  along  the  z  axis.  ’I  he  curves  have  liven  dis- 
piaeed  to  the  U  it  and  the  leading  edges  made  lo  coincide  at  ' ' 
height  lunu  (erg  emJ),/*  lor  I  £  I  .  the  abscissa  For  this  Inter¬ 
section  of  the  curves  has  been  lulieied  42.  a,  the  location  ol 
this  (Hilnl  at  t  0. 

linear  systems.  This  argument  is  presented  in  detail  in 
another  paper  and  another  kind  of  criterion,  "utility", 

Is  proposed.  A  difference  scheme  is  "useful"  for  pre¬ 
scribed  (  and  S'  if  the  contpuicr  solution  of  the  scheme 
produces  results  for  S'  time  steps  which  differ  In  mag¬ 
nitude  front  the  correct  solution  of  the  original  differen¬ 
tial  equations  by  an  amount  less  than  «.  Thus,  utility 
analysis  is  totally  different  in  philosophy  titan  stability 
anti,  Incidentally,  is  far  easier  to  apply.  A  number  of 
theorems  have  been  developed  which  allow  rapid  and 
powerful  assessment  of  the  utility  regions  of  differencing 
schemes  of  all  types.  Because  tiiese  matters  would  only 
extend  the  length  of  Iho  present  paper  and  because  it  Is 
believed  that  utility  procedures  wiii  be  worthwhile  In  a 
braider  context  than  the  problem  at  hand,  a  separate 
paper  Is  devoted  to  the  subject  and  only  a  brief  sketch  of 
the  analysis  will  be  given  here. 

It  turns  out  that  a  utility  region  for  an  explicit  differ¬ 
encing  scheme  can  he  determined  from  a  relation  of  the 
form 

A/<1(46)(AtA')]-,1  (6) 

where  A t  Is  the  time  step  size,  <V  is  the  number  of  time 
steps  to  be  made,  f  ~10'J,  and  Z  Is  a  function  of  the  spa¬ 
tial  mesh  sizes  and  Is  determined  from  the  differencing 
seneme.  For  the  seven-point  explicit  difference  scheme 
used  to  represent  the  laser-fnid  system,  6  Is  Irlvialiy 
determined  to  be1 

Z*c\3/kL(r,j&xf  +  2/Az],  (7) 

so  that  the  utility  restriction  prcdlcled  is 


-y-  (r  »  200cm) 
0  ° 


FIG.  U.  lladlal  profile  of  Ifl  la  shown  for  t-0  and  for  t-  Id"4 
see  for  slices  taken  through  the  on-axis  maximum  z ^  in  the  z 
profile. 
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K1G.  4,  Details  of  the  radial  profile  are  shown  for  various 
times.  In  all  four  eases  the  radial  slice  through  the  on-axls 
maximum  of  the  z  profile  is  exhibited. 


As  a  start  ng  point  for  this  discussion,  the  parameter 
values  used  in  the  actual  calculation  will  be  listed.  The 
electric  field  at  f  =  0  was  taken  to  be  of  the  form 

£,  =  Fexpf-^r/Vo)2]  exp{-  4[(2  -  zj/zj] , 

£t-°f  (9) 

where 

r0«fuil  (1/e)  width  of  Et(r,z  =  z„  1  =  0) 

■/2 1  full  0/e)  width  of  IL  at  z=z0O,  f  =  0); 

4  *  [on-axls  Intensity  (in  erg/cm'sec)  at  r  =  O,z  =  zp0, 

1  =  0,  time  averaged  over  several  optical  periods] 

■Ifco .WcF': 

*0  =  full  (1/e)  width  ol  £,(7  —  0,2,  f  =  0). 

*,0=  location  of  the  peak  at  f  =  0 ; 

F2  21  l  .  16/1 _ 64  U 

*  if (2»)‘  '\,r^z„  ^ 

=  peak  value  of  the  electric  Held,  squared; 

PL  =  total  power  (In  erg/sec)  of  the  pulse  at  t  =  t^, 

1  =  0,  time  averaged  over  several  optical  peri¬ 
ods; 

U=  total  energy  In  the  pulse  at  '  =  0,  time  averaged 
over  several  optical  periods. 

The  values  taken  for  these  quantities  were 
r0-  200  cm , 

20  =  9  xlO5  cm  =  9  km  , 

2^,  =  13.  5x10’  cm  =13  km  , 

£=4.7x10’  (erg/cm3),/,  ,  (11) 


4  =3.  3  x  lo”  erg/(cm2  sec) , 

PL  =  5.  2X102'  erg/sec  , 

{/=  9.  8  x  10"  erg, 

and  the  air  was  taken  to  he  Initially  in  Its  unperturbed 
state  at  1  atm  pressure  and  at  10 *C.  The  reason  for 
these  astronomical  powers  will  be  explained  later. 

The  spatial  grid  was  composed  of  80x15=  1280  mesh 
l»ints.  The  z  axis  was  evenly  divided  into  80  steps  of 
size  A2  =  0. 45xl05  cm  =0. 45  km  beginning  at  2  =  0  and 
extending  to  z  =  35.  6  x  10’  cm  =  35.  6  km.  Thus  the  peak 
of  £,  was  initially  located  at  the  30th  mesh  point  on  the 
2  axis  and  its  1/e  width  extended  from  the  20th  to  the 
40th  mesh  point.  The  radial  variable  x  has  the  range 
Oix£l  and  this  range  was  evenly  divided  into  22  steps 
of  size  Ax  =  ^j,  but  oqly  the  16  sites  closest  to  the  z 
axis  were  used.  The  more  distant  sites  correspond  to 
radial  distances  greater  than  five  beam  half-widths.  The 
first  step  away  from  the  z  axis  corresponds  to  the  radi¬ 
al  distance  Ar=  200  cm/21  =9.  5  cm  *^x(radial  half- 
width).  The  time  step  size  was  taken  to  be  A/^IO"1  sec 
and  100  steps  were  made  so  that  the  time  interval  0  <,t 
<  10'5  sec  was  traversed. 

Taking  these  grid  sizes  and  time  steps  and  subtitutlng 
into  the  utility  condition  (8),  one  gets 


-f-  f  r„  »  2  00cm )  » 

f0  0 

FIG.  5,  An  off-axis  maximum  in  the  radial  profile  is  shown 
at  /-  10*5  sec.  The  slice  shown  exhibits  the  radial  profile  at 
2  31,  whereas  die  principal  peak  of  the  pulse  is  on  axis  at 
2  =  35.3.  The  slice  at  2  =  31  contains  the  greatest  off-axis  effeet 
and,  therefore,  locates  the  two  secondary  peaks  whieh  have 
developed  In  the  pulse.  These  secondary  peaks  are  aiso  indi¬ 
cated  on  Fig,  9. 
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FIG.  G.  Phase  Information  Is  presented  lay  showing  I  C(  I  as  a 
functlnn  of  z  at  I  Ux  Id*’'  sec.  l-c»r  comparison,  llie  dashed 
curve  shows  I  £1  .  K  (illations  14  iuk!  (14  of  the  .ext  define  |£,| 
and  I £ I  . 
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cA/~  3000-40800  5030 

=  l2/kLr§(M?  +  8/ixz.  (12) 

Thus,  for  the  chosen  step  size  Az,  the  value  used  for 
A/  would  violate  the  ccndltion,  were  it  to  be  doubled.  Of 
course,  the  condition  (8)  Is  only  approximate,  but,  as 
mentioned  above,  pood  solutions  could  not  be  obtained 
for  A/~10  l  for  Az  =  0.45  '  i,  Clearly  the  value  fa  used 
for  Ax  does  not  saturate  the  Ax  piece  of  (12)  and  one 
could  probably  use  Ax  as  small  as  ^7.  Such  a  small  step 
size  for  Ax  would,  however,  require  three  times  as 
many  spatial  mesh  points  and  would  exceed  the  storage 
capacity  of  the  computer  which  was  used. 

The  desire  is  to  use  as  large  a  value  of  Af  as  one  can. 

In  this  regard,  the  Ax  piece  of  (12)  Is  generous  and 
would  permit  A/-10’".  The  z  step  size  would  have  to  be 
Increased  to  Az~5  km  to  allow  this,  however.  Such  a 
large  step  size  would  be  larger  than  the  4.  5  km  half- 
width  of  the  pulse  selected,  so  that  no  details  of  defor¬ 
mation  of  the  pulse  could  be  observed. 

If  the  beam  Is  made  narrower  in  radial  extent,  the 
Fresnel  length  decreases  and  diffraction  effects  become 
Important.  The  Fresne.  length  is  310  km  for  ro=200  cm, 
so  that  one  would  become  Involved  with  far-fieid  effects 
If  the  beam  radius  were  decreased  by  more  than  a  factor 
of  5.  Making  the  pulse  longer  in  tiie  z  direction  expands 
the  time  scale  over  wiiich  interesting  effects  may  be 
studied.  If,  on  the  other  hand,  the  pulse  is  shortened  in 


the  z  direction,  then  one  must  shift  to  smaiier  values  of 
Az  In  order  to  h"  able  to  follow  details  of  the  develop¬ 
ment  of  the  pulse.  Shifting  to  smaiier  Az  requires,  be¬ 
cause  of  (12),  that  one  use  smaller  values  of  A /  The 
net  effect  is  that  no  profit  is  derived  from  using  shorter 
pulses,  because  they  ran  be  followed  only  for  corre¬ 
spondingly  shorter  t/rhes. 

One  aspect  of  the  parameter  regime  has  not  yet  been 
discussed:  the  range  of  power  for  the  beam.  Since  the 
power  4  does  not  appear  in  the  utility  criterion,  Its  role 
must  be  determined  by  experimentation  with  the  comput- 
ei  program.  Ve\-f  small  powers  are  not  interesting  be¬ 
cause  there  is  very  liltle  interaction  with  the  fluid.  In 
order  to  see  instabilities  and  nonlinear  effects  during 
short  times,  one  would  wish  to  consider  beams  with 
large  power  densities.  The  extremely  large  values 
shown  in  (11)  produce  interesting  effects,  In  a  time  in¬ 
terval  of  10'5  sec.  Such  beams  cannot  be  followed  for 
more  than  about  100  time  steps,  however,  because  the 
various  dependent  variables  begin  to  develop  large  cur¬ 
vatures  and  vary  on  a  scale  smaller  than  the  mesh  sizes. 
Thus  If  one  wishes  to  follow  the  development  for  a  long 
period  of  time,  the  mesh  sizes  must  be  decreased  and 
eventually  the  time  step  will  have  to  be  smaiier,  and 
then  many  more  time  steps  will  be  required.  In  this  re¬ 
gard,  one  must  keep  in  mind  that  if  the  mesh  size  Is  de¬ 
creased,  while  the  init.-il  puisc  size  is  not  decreased, 
then  more  mesh  points  will  be  required  and  the  storage 
capacity  of  the  computer  also  becomes  a  limiting  factor. 


FIG.  7.  At  t» 8x  10*®  see,  I £(  I  and  l£i  are  shown  on  axis  as 
functions  of  z.  Two  nodes  have  developed  and  £,  Is  negative  In 
the  region  of  the  |jower  peak.  The  sign  ol  £(  in  the  various  re¬ 
gions  is  Indicated  on  the  figure.  The  nodes  are  also  shown  in 
Hg.  9. 
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I'lCi.  it.  At  /  1  ()■'■'  sec,  IKjl  and  I  A' I  arc  shown  on  axis  as 
functions  of  z.  The  sign  of  f,  is  indicated  in  the  various  re¬ 
gions.  Thors'  art'  now  four  nodes.  The  nodal  carver  arc  phase 
fronts  and  arc  shown  in  detail  in  I  Ik.  9. 

The  final  remaining  option  is  to  increase  the  |>ower  in 
tholio.ini  oven  more.  The  not  oifocl  is  t tint  tho  large 
curviituros  develop  fastor  and  tho  development  can  bo 
foiiowed  only  for  shorter  periods  of  tinio. 

One  finai  remark  about  numerical  solution  of  tho  iasor- 
fiuid  equations  will  bo  made  before  discussing  tho  re¬ 
sults  of  tho  computer  calculation.  St rong  growth,  insta¬ 


bilities,  and  nonlinear  effects,  can  often  not  be  followed 
because  of  the  mesh  sizes  employed  If  these  ‘troop 
oscillations  or  secular  growths  are  pone  rated  by  tiny 
rapidly  dunging  terms,  that  is.  if  the  instabilities  arise 
due  to  ripple  effects  which  become  strongly  enhanced, 
then  a  crude  mesh  size  can  smooth  these  effects  out  and, 
thereby,  prohibit  the  occurrence  of  the  strongiy  growing 
phenomena  l>y  removing  their  source  Very  strong  in¬ 
stabilities  were  found  in  Hof,  1  for  the  lineal  ized  i.iser- 
fluid  equations.  The  strongest  of  these  instabilities  are 
generated  by  very  short  wavelength  ripple.  The  mesh 
size  employed  in  the  present  calculation  wiil  begin  to 
wash  out  rip, pie  altout  an  order  of  magnitude  iarger  in 
wavelength  than  the  ripple  which  is  most  strongly  ampli¬ 
fied  in  the  linearized  analysis.  Thus,  one  must  bear  in 
mind  that  some  physical  sources  of  pulse  dbtortion  wiii 
be  excised  by  the  me.sh  seieeted. 

Accepting  the  many  restriction.!  noted  above,  we  have 
examined  the  propagation  of  a  200-rm  by  9-km  puise 
with  10IT  erg  for  10'"‘  see.  Tho  pulse  moves  3  km  during 
this  time  and  it  is  possible  In  observe  the  onset  of  the 
iTser-fiuid  interaction  In  some  detaii. 

111.  RESULTS  OF  THE  COMPUTATION  FOR 
THE  ELECTRIC  FIELD 

The  results  of  the  calculation  are  presented  in  Figs. 
1—17.  Tlie  eieetric  ficid  is  conveniently  considered  in 
terms  oi  the  quantify 

|£|  .(£*4  £*)•/*,  (13) 

where  /■',  ”nd  A,  are  the  slowly  varying  electric  ampli¬ 
tudes  defined  in  (2).  The  instantaneous  eieciric  fieid  is 
thus  given  by 

£>|E|cos(.et,/-*i>z4  M,  04) 

where  the  phase  bB  is  given  by 

As  shewn  in  (9),  at  /  — 0,  l  2  it  taken  to  be  zero  and, 
consequently,  br  Is  zero  inltiaiiy.  Thus 
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I’K*.  9.  Various  properties  of  the  pulse 
are  .-dmwn  in  tiie  rz  plane.  The  location  of 
the  peak  in  the  z  profile  is  shown  as  a 
function  of  r  at  /  =  ()  and  at  /  10-5  see.  The 
phase  fronts  with  E,  -0  are  shown  at  I 
t  10*5  see.  Tiie  open  circles  locate  tin  z- 
profile  nodes  of  f,  at  /  -HylO"*  sec.  The 
small  squares  locate  the  secondary  max¬ 
ima  of  the  pulse  at  t-  It)'11  sec. 
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HC..  lo.  On-axis  loinpe  ratine  increment  T  -  T„  la  shown  as  a 
fund  Ion  of  z  at  I  lit'1  see. 


so  that  after  10  sec,  the  p»ak  has  lost  about  J  kr.i  with 
respect  to  zf.  Note  that  the  exponential  damping  factor 
shown  In  (2)  is  not  included  u\  the  quantity  1 1  |  appearing 
in  the  graphs.  For  air,  tills  factor  Is  iarger  than  0.95, 
even  at  /=  10  sec.  Other  than  this  effect,  very  little 
energy  is  lost  from  the  beam  due  to  heating  of  tiie  fluid, 
so  the  distortion  effects  shown  in  Fig.  1  are  rather 
minor  and  are  noticeable  only  near  the  peak  of  the  pulse, 
bxtra  detail  of  this  peak  distortion  1.,  shown  in  Fig.  2. 
For  purposes  of  tills  display,  tiie  leading  edges  have 
been  placed  together  so  that  the  curves  intersect  at  I /•  | 

=  10*  (crg/cm1)' /2  and  tiie  corresponding  abscissa  lias 
been  labeled  42.  5,  the  location  of  this  point  at  /  0.  The 
retrograde  peak  motion  and  the  corresponding  loss  of 
fore  and  aft  symmetry  in  the  vicinity  of  the  peak  are 
plainly  seen. 

The  radial  beam  profile  is  exhibited  in  Figs.  3  and  4. 

Tiie  radial  slice  shown  at  each  value  of  the  time  is  taken 
through  tiie  position  of  the  maximum  z.  in  the  z  profile. 

In  Figs.  3  the  entire  beam  profile  is  shown  for  the  ini¬ 
tial  and  final  times  only.  Comparison  of  the  curves  re¬ 
veals  a  small  on-axis  increase  extending  out  to  the 
beam  half-width  (half  of  the  full  1  /c  width)  at  !.ib=  100 
cm.  Tiie  effect  amounts  to  a  47'j?  increase  in  the  on-axis 
intensity.  iX*taiis'  and  intermediate  states  are  given  in 
Fig.  4. 

Although  it  may  appear  from  Figs.  1-4  that  the  energy 
In  the  pulse  is  not  conserved  properly  and  that  the  beam 
is  gaining  energy,  such  is  not  the  case.  The  total  energy 
in  the  pulse  is  in  fact  constant  to  within  5r?  throughout 


/•'I  =  A',  at  /  —  0 , 


(15) 


and  E,  is  described  by  Kqs.  (9)-(l  1 )  initially.  This  ini¬ 
tial  pulse  shape  in  exhibited  in  Figs.  1(a)  and  3. 

In  Fig.  1  the  z  profile  of  tiie  pulse  is  shown  at  the  initial 
time,  at  10'5  sec,  and  at  two  intermediate  times.  For 
1*0,  the  pulses  arc  not  absolutely  symmetric  atiout 
their  peaks.  In  order  to  exhibit  this  asymmetry,  tiie 
curves  arc  plotted  as  a  function  of  I  z  -  ze] ,  where  z  is 
the  center  of  ‘.lie  pulse.  Tills  device  allows  direct  com¬ 
parison  of  the  leading  and  trailing  edges  of  tiie  puises. 
The  center  zc  is  defined  to  lie  the  point  equidistant  from 
the  leading  and  trailing  edges  at  I A I  =  1  (crg/cm5)1 /J. 
Tlies  ;  values  are  (/  in  sec) 


f  =  0  f=r  G  x  10'*  f  =  8xi0'*  f  10*’ 

z(  =  30  ze*34  zc  *35. 3  zc*36,6, 


(16) 


where  for  convenience,  distances  along  the  z  axis  will 
be  given  in  units  of  grid  size:  Az=.  0.  45  x  10’  cm  =  0,  45 
km.  One  notes,  therefore,  from  (16)  that  this  puisc  cen¬ 
ter  propagates  at  the  velocity  rc  *2.  97x  10' 1  cni/sec,  the 
velocity  of  light.  The  pulse  peaks,  iiowevcr,  are  ob¬ 
served  to  drift  backward  with  respect  to  zc  (t  In  sec): 


/  =  0 


z^  *30 


z,*34 


6x10'*  /  =  8 x  10'6  /=  10'5 
*35. 3  , 


(17) 


z,* 35.  2 


I*.*  ?OOcm»  - 


!??’  Ua,Jial  Profile  of  the  temperature  Increment  is  shown 
“l. 1  1,“  sec  for  l,MS  s,ieL'  through  tii-  maximum  of  the  z  pro- 
°*  *iKL\imum  *  S'A  as  may  be  seen  in  i  ig.  lo. 
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1'IG.  12.  i  component  r,  of  (he  fluid  velocity  Is  shown  on  axis 
as  n  function  of  z  at  I  10"s  sec,  A  double  Ion  plot  Is  used  w  hich 
omits  values  of  r,  between  10*5  and  -  10*5  cm/sec. 


the  computer  calculation.  The  shape  of  the  pulse  Is 
merely  becoming  slightly  more  complicated.  Although 
the  radial  peak  is  on  axis  In  the  slice  through  the  peak 
in  the  z  profile,  this  Is  not  the  case  for  slices  taken  be¬ 
hind  zf.  For  example,  at  / trlO"’  sec,  the  principal  peak 
is  at  zf  *35.  3.  As  one  moves  away  from  this  peak  to¬ 
ward  the  trailing  edge,  the  radial  peak  moves  off  axis 
giving  a  maximal  effect  near  4  -31.  The  principal  peak 
of  the  pulse  is,  however,  always  on  axis.  The  radial 
profile  at  4  =  31  is  shown  in  Fig.  5.  This  is  clearly  only 
a  small  detail  at  /  — 10"’  sec.  The  position  of  this  off- 
axis  secondary  peak  is  also  located  in  Fig.  9  and 
marked  with  tiny  squares. 

In  order  to  follow  the  development  of  the  phase  of  the 
electric  field,  the  quantity  I  A',  I  is  plotted  in  Figs.  6—8, 
for  1*0.  Of  course,  at  t  —  0,  I  A,  I  —  |  A'  I  and  the  phase  bE 
Is  zero  In  these  three  figures  the  graph  of  I  A I  Is 
marked  with  dotted  lines  for  comparison.  The  corre¬ 
sponding  value  of  I  A, I  can  be  deduced  from  these  fig¬ 
ures,  using  (13).  These  figures  show  far  more  dramatic 
effects  than  the  curves  discussed  above.  At  f=6-><10'" 
sec,  the  phase  is  still  nearly  zero  and  A,  is  positive 
everywhere.  At  ^  =  8 >■  10-"  sec,  however,  A,  has  changed 
sign  over  a  3-km  region  extending  from  slightly  in  front 
of  the  peak  of  1 A I  toward  the  trailing  edge  of  the  beam. 
This  is  clear-cut  evidence  of  the  onset  of  laser-fluid 
interaction  In  the  trailing  edge  of  the  beam.  It  Is  clear 
that  the  front  ol  the  pulse  and  the  distant  tail  are,  as 
vet,  unaffected  by  this  interaction,  it  is  interesting  that 
the  unperturbed  part  of  the  leading  edge  does  not  reach 
as  far  back  as  the  principal  peak.  Thus,  the  peak  al¬ 
ready  feels  the  effects  of  the  Interaction  to  some  degree 


and  comparison  of  (16)  and  (17)  reveals  that  tiie  peak 
will  now  begin  to  lose  ground  with  respect  to  the  center 
of  the  pulse.  This  effect  has  already  been  noted  in  the 
graphs  of  I  £1 . 

Since  there  are  now  places  where  I  A,  I  Is  zero,  it  is 
clear  that  the  phase  goes  to  Jit  at  these  sites.  The  am¬ 
plitude  A,  responds  strongly  at  those  places  where  A, 

=  0,  fulfilling  the  obligation  to  conserve  power.  One 
notes  that  the  graph  of  I  A'!  remains  very  smooth,  giving 
no  indication  that  the  phase  is  varying  rapidly.  Figure  8 
shows  the  later  development  of  the  region  in  which  A, 
changed  sign.  The  region  in  which  A,  <0  is  now  7  km 
long,  nearly  as  large  as  the  1/c  width  of  IAI.  This  re¬ 
gion  has  advanced  now  to  point  well  in  advance  of  the 
principal  peak  and  extends  back  far  into  the  tall.  It  ap¬ 
pears  that  this  node  is  propagating  forward  at  nearly 
four  times  the  speed  of  ligiit.  Furthermore,  there  has 
been  another  sign  revere  ti  of  A,  slightly  behind  the 
peak.  It  is  this  kind  of  oscillatory  behavior  in  A,,  with 
large  variations  on  the  scale  of  the  chosen  mesh  size, 
that  brings  a  halt  to  further  observation  of  the  beam  de¬ 
velopment  by  this  method. 

On  the  rz  plane  shown  in  Fig.  9,  the  constant  piu  se 
curve  A,  =  0  is  shown  in  detail  at  /.=  10‘*  sec.  Also 
marked,  with  small  open  circles,  is  the  on-axis  extent 
of  the  similar  curve  at  f  =  8*10‘"  sec,  encountered  in 
Fig.  7.  Also  indicated  on  the  same  figure  is  the  locus  of 
maxima  in  the  z  profile  of  the  pulse  for  off-axis  slices, 
both  at  /  =  0  and  /  =  !0*'  sec.  At  /  =  0,  the  pulse  is  de¬ 
scribed  by  (9)  and  clearly  the  off-axis  slices  all  have 
maxima  in  z  positioned  at  4^  =  30.  As  the  pulse  propa- 
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FIG.  13.  Radial  component  rr  of  the  fluid  velocity  Is  shown  as 
a  function  of  r  at  t~  10"5  see.  The  slice  Is  taken  at  *  =  32,  the 
location  of  the  “center  of  velocity”  shown  In  Fig.  12. 
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£=^(r,z-r/)exp{4<i>tf-(v7w1/r)z]}  +  c.c  ,  (18) 

where 

^0»  ifexp[-4(rA0)2]exp{-4[(z-c/)/z0P  (19) 

and 

v7*l+[(fo-l)/2po]ft,  (20) 

wl  ere  F  is  a  slowly  varying  amplitude,  p,  =  p-po  is  the 
local  density  excess,  and  for  ttie  present  considerations 
ip o  has  been  put  equal  to  zero.  Combining  Eqs.  (18)— 
(20),  the  electric  Held  can  be  put  in  the  form  shown  in 
(2)  with 

£- ■[  ■ f  cos(‘“  ‘Sir')] i  sT]  “p[■''(£^^),]  • 

(21) 

These  expressions  agree  with  (9)  at  I  =  0  and  offer  a  way 
estimate  the  behavior  of  £,  at  subsequent  times.  On  the 
basis  of  (21)  the  nodes  of  £,  might  be  expected  to  be  de¬ 
termined  by 

cos^zUfo-D^Polp.JsO.  (22) 

Actually,  this  expression  should  be  modified  slightly  if 
one  wishes  to  attempt  to  get  quantitative  agreement  with 
the  computer  solution.  It  is  clear  that  (21)  requires  Et 
to  vanish  at  z  =  0  at  all  times.  This  is  not  the  same 
boundary  condition  which  was  u  .ed  in  the  computer  solu¬ 
tion.  Actually  one  should  use 


gates,  however,  the  peak  moves  slower  than  the  velocity 
of  light,  as  previously  noted.  The  off-axis  portions  of 
the  pulse,  however,  have  much  smaller  intensity  and, 
consequently,  interact  very  little  with  the  fluid.  These 
portions  of  the  pulse  will  suffer  no  delay,  and  move 
steadily  ahead  of  the  principal  peak.  One  notes  that  at 
two-thirds  of  the  radial  half-width  of  the  beam,  the  de¬ 
lay  has  disappeared  almost  completely.  As  mentioned 
above,  the  small  squares  locate  the  secondary  peaks 
present  at  l=  10"s  sec. 

IV.  analytical  procedure  developed 

TO  EVALUATE  AND  EXTEND  THE  CALCULATION 

The  phase-front  information  presented  in  Figs.  6—9  is 
an  interesting  feature  of  Ihe  results  of  the  computer  so¬ 
lution.  Since  |£[l  turns  out  to  be  a  rapidly  varying  func¬ 
tion  of  time,  it  is  of  interest  to  attempt  to  understand 
the  mechanism  responsible  for  the  behavior  of  £,.  In 
order  to  understand  this  behavior,  one  must  realize  that 
the  phase  depends  on  the  state  of  the  fluid.  The  state  of 
the  fluid  given  by  the  computer  calculation  Is  shown  in 
Figs.  10—16.  Before  these  figures  are  discussed,  how¬ 
ever,  it  is  convenient  to  examine  certain  analytic  esti¬ 
mates  for  the  fluid  variables.  Such  estimates  will  allow 
insight  into  the  behavior  of  £,  and,  later,  will  facilitate 
the  discussion  of  the  computer  results  for  the  fluid 
variables. 

In  order  to  describe  the  behavior  of  £,,  it  is  useful  to 
write  the  electric  field  in  the  form 


£,  ~CO'ilp 

with 


IO 


FIG.  15.  Radial  density  distribution  Is  exhibited  as  a  function 
of  r  at  I*  10"*  sec.  The  slice  Is  taken  at  i  -  Si,  the  location  of 
the  density  minimum  detailed  in  Fig.  14.  A  double  log  plot  is 
used  which  omits  values  lietwcen  It)"’  and  -10*’  g/cni  . 
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too 


KiG.  It;.  On  a  double  log  plot,  the  vntlous  laser  and  fluid  vnrl- 
alilos  am  simultaneously  plotted  versus  z  at  /  10*5  sec  so  that 
the  spatial  location  of  the  various  pulses  can  he  visualized. 


•i~hL  (23) 

The  integrand  in  (23)  should  be  evaluated  at  i'  and  /', 
where 


cW -l)  =  z' -z. 


It  is  clear  from  (21 )  and  (23)  that  the  mechanism  re¬ 
sponsible  for  the  behavior  of  £',  is  easily  exhibited.  To 
actually  follow  the  behavior  of  however,  it  is  cieariy 
necessary  to  determine  the  state  of  the  iluid.  In  particu¬ 
lar,  the  density  excess  p,  must  be  obtained  as  a  function 
of  time  and  position.  In  order  to  analytically  describe 
the  fluid  for  the  time  interval  and  parameter  ranges  of 
the  computer  solution,  the  laser-fluid  equations  may  be 
simplified  to 


A> 


3/  (3  3/  t„ 


(24) 


and 


jp-'u’fov’T, ,  (25) 

where  7|  is  the  local  temperature  excess,  T-T0.  One 
must  recall  that  the  intensity,  cP,  appearing  in  (24)  is 
a  function  of  time  and  position.  Integrating  (24)  from 
zero  to  /  and  combining  the  resulting  equation  with  (25) 
to  eliminate  Tx,  one  can  obtain 


It  is  straigiitforwai  d  to  integrate  this  expression  and, 
although  the  details  wiii  not  be  given  here,  it  is  clear 
that  p,  will  have  the  form 


Pi 


[(^-l)exp(-^l)]/-U,/),  (27) 


so  that  p*"  p„  on  axis  and  there  is  an  off-axis  maximum 
in  the  density.  In  other  words,  there  wiil  be  a  piie  up  of 
the  fluid  at  a  distance  r*{ra  from  the  axis. 


The  temperature  distribution  can  easiiy  be  obtained  bv 
integrating  (24)  from  zero  to  /: 


71* 


A  negligible  term  involving  p,  can  be  evaluated  using 
(27)  and  has  been  dropped  to  obtain  (28). 

Combining  (27)  will  (23)  the  behavior  of  Ex  can  be  visu¬ 
alized  and  studied  analytically.  Closed  contours  «■  ich  as 
those  shown  in  Fig.  9  are  predicted  and  other  quaiit  Jive 
features  are  correct.  A  detailed  comparison  of  this  an¬ 
alytical  procedure  and  the  computer  resuit  is  in  pro¬ 
gress  and  it  is  now  clear  that  striking  quantitative 
agreement  Is  obtained  This  succ  ss  is  of  great  Interest 
because  the  analytic  procedure,  i.niike  the  present  com¬ 
puter  solution,  is  not  limited  to  10'’  sec.  For  the  pres¬ 
ent, 'however,  only  the  computer  solution  is  presented 


V.  RESULTS  OF  THE  COMPUTATION  FOR 
THE  STATE  OF  THE  FLUID 


Analysis  of  the  computer  result  for  the  state  of  the  fiuid 
at  / =  10’’  sec  wiii  now  be  considered.  For  this  discus- 
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PIG.  17,  Location  and  full  1  fc  widths  of  the  various  laser  and 
fluid  pulses  arc*  shown  vs  z.  The  |x*ak- to- valley  distance  is 
shown  for  i'#.  The  Initial  and  finai  location.'  of  the  laser  power 
peak  are  also  indicated. 
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sion,  it  is  useful  to  keep  severai  ciiaracteristic  dis¬ 
tances  in  mind,  Since  the  state  of  the  fluid  is  governed 
by  the  tntensity  profile 

IL{r,z,l)*':((0')'/1c\l{2,  (29) 

rather  than  by  /■  directly,  the  following  initial  param¬ 
eters  are  relevant:  The  full  l/o  width  of  the  z  profile  of 
the-  intensity  ts 

Zo/.'jj  =  1 4. 1  in  units  of  d£  « 0  45  km.  (30) 

At  £  =  £,±  \x0-z,,±  5, 


so  that  ?lL/?z  has  extrema  separated  by  10  units.  The 
full  -J-  width  of  the  radial  intensity  profile  is 

»o/v2  =0.70Vo=141  cm.  (32v 

At  r  =  1  r0  =  50  cm, 


(33) 


so  that  lL/br  has  a  maximum.  The  maxima  on  opjjosite 
sides  of  the  axis  fdi.mieter  of  the  density  doughnut)  are 
separated  by  100  cm, 

* 

The  local  temperature  excess  T-Tn,  where  7(,=  10*C, 
is  shown  in  Fig  10  as  a  function ot  z  at  /-10"1  sec. 

This  curve  is  in  complete  quantitative  agreement  with 
the  analytical  result  shown  m  (2m).  The  hottest  plaee  in 
tin  beam  lies  on  the  axis  at  £=33.3.  Since  the  intensity 
peak  is  at  zf  *35.3,  it  is  clear  that  the  thermal  peak  is 
lagging  behind  the  intensity  peak.  Since  £r  =  36.  G  at  this 
time,  it  is  clear  that  the  thermal  peal;  is  almost  exactly 
midway  between  the  initial  and  final  pulse  centers.  Thus 
one  finds,  as  expected,  that  the  thermal  peak  propagates 
at  velocity  |.r  tor  small  times.  This  and  oilier  proper¬ 
ties  of  the  thermal  profile  are  readily  understood  on  the 
basis  of  the  lidiowing  considerations.  Tin.  temperature 
responds  to  the  heat  deposited  in  the  medium,  so  that 
from  (28), 


nt)-ra*f‘  (~-iLu,z\i^ oyif) .  04 > 


Due  to  the  symmetry  of  lL  and  the  fact  that  it  is 
Gaussian,  it  follows  from  (34 )  that  T(z)  should  reach  a 
maximum  midway  between  z^  and  zp  for  small  times. 
Furthermore,  the  grapli  of  7(z)~  T„  sli.  utd  lie  symmet¬ 
ric  about  its  maximum.  Both  of  these  f  'atures  are  evi¬ 
dent  in  Fig.  10.  Since  the  peak  moves  much  less  than 
Us  half-width  in  10'  sec,  the  integrand  in  (34)  is  essen¬ 
tially  constant.  Kvaiuating  lL(y,z’,l)  at  the  mid|ioint  of 
the  interval,  one  obtains 

T(z)-  To'HZ'-z^/cHn/faCjl^r,  z  -  (( zp-  z^),  <  =  0). 

(35) 

From  (35)  one  concludes  that  the  width  of  the  thermal 
distribution  should  equal  the  width  of  the  intensity  dis¬ 
tribution.  Indeed,  one  sees  lit  Fig.  10  that  the  tempera¬ 
ture  distribution  has  a  width  14,1,  which  Is  to  lie  com¬ 
pared  with  (30). 


The  corresponding  radial  temperature  distribution  at 
£-33,  the  position  of  the  maximum  at  /=10*5  sec,  ts 
given  in  Fig,  11.  Tills  curve,  also,  ts  In  complete  quan¬ 
titative  agreement  witli  the  analytical  result  shorn  in 
(28).  The  temperature  lias  reached  a  maximum  of  more 
than  1 000  °K  on  axis  and  the  full  width  of  the  distribution 
is  lound  to  be  0,  67>-0  =  134  cm,  5ra  narrower  than  the 
initial  radial  intensity  width.  One  mig’ t  have  expected 
the  temperature  d'sM  iimtion  to  lie  broader  than  the  in¬ 
tensity  profile  because  the  targe  radial  velocity  of  the 
fluid  should  carry  some  of  the  deposited  energy  away. 

In  fact,  this  effect  may  jiossibly  lie  observed  in  the  fol¬ 
lowing  way.  One  might  compare  the  temperature  distri¬ 
bution  not  to  the  original  Gaussian  intensity  profile,  but 
rather  to  the  1/c  widtli  of  the  radial  intensity  profile  at 
/=  10"’’  sec.  This  final  intensity  profile  has  a  width  of 
128  cm  or  5  less  than  the  temperature  width.  The 
average  of  the  two  intensity  widths  is  134.  5  cm,  aimost 
exactly  the  observed  temperature  width.  Tills  average 
may  ire  the  best  measure,  bccuase  the  thermal  peak  is 
midway  between  (lie  ini''  i  and  final  intensity  peaks. 

Tiie  z  component  ot  the  fluid  velocity  Is  shown  in  Fig. 

1 2  as  a  function  of  z.  A  double  log  plot  is  used  which 
omits  values  of  i/  between  10'1  and  -10'1  cm/see.  This 
kind  of  plot  allows  negative  values  of  rI  to  be  plotted 
beiow  the  “axis”.  The  zero  of  the  velocity  distribution 
occurs  around  £  =  32.  2,  so  the  “center  of  velocity”  lags 
slightly  behind  the  thermal  maximum.  From  the  differ¬ 
ential  equation  for  i/t  one  might  expect  to  find 


’  (36) 

so  that  the  peaks  in  Fig.  12  would  be  separated  by  10 
units  according  to  (31).  Indeed,  the  peaks  are  found  to 
be  separated  by  10.  2  units.  Furthermore,  since  the 
temperature  curve  is  symmetric  about  its  maximum, 
(36)  would  suggest  that  should  lie  antisymmetric  about 
its  zero.  Tills  effect  is  correctly  observed  in  Fig.  12. 
The  velocity  distribution  is  delayed  with  respect  to  the 
temperature  distribution,  but  this  symmetry  property  is 
unaffected. 

At  the  same  value  of  z,  corresponding  to  the  center  of 
velocity  of  the  z  component,  the  radial  component  vr  is 
plotted  in  Fig.  13.  This  value  of  z  corresponds  also  to 
the  largest  radial  velocities,  so  that  £  =  32  might  also 
be  termed  the  site  of  greatest  kinetic  energy  in  the  fiuid. 
The  radial  velocity  maximum  is  4600  times  larger  than 
the  maximum  axial  velocity.  This  effect  arises  because 
of  the  great  disparity  in  the  intensity  widths  in  the  two 
directions.  As  a  matter  of  fact, 


-|^i=4625  and 


9  xlO5 
200 


=  4500. 


The  curve  of  i»f  Is  forced  to  go  to  zero,  as  r  goes  to 
zero,  by  the  boundary  conditions  shown  in  (5).  One 
notes,  however,  that  the  peak  is  located  at  the  distance 
0.  25r0  from  the  axis,  exactly  the  location  of  the  maxi¬ 
mum  of  dlL/br  shown  in  (33).  Thus  one  finds 


<t/,  ?T 

V  ‘ - *“ - 

r  fV  ?r 


(37) 
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as  would  bo  expected  from  the  diffcrentail  equation  for 

'r 

Flringing  up  the  reir  in  tlie  sequence  of  effects  is  tiie 
density  minimum  it  **31.5  Tiie  density  decrement 
-  (p  -  f\„)  is  shown  in  Fit;.  14  Using  (25),  one  mighl 
expect 


?  P  ST*  ’ 


(38) 


so  that  p  anti  T  will  have  the  same  *  dependence.'1  There 
should,  however,  lie  a  double  time  delav,  since  two  time 
integrations  are  indicated  In  (38).  The  width  of  tiie  den¬ 
sity  decrement  Is  found  to  lie  14.  5.  about  3fJf  wider  than 
the  thermal  and  intensity  widths  The  maximum  frac¬ 
tional  decrement  (p„-  p)/p0  ~  4  *  10*’. 

The  radial  density  distribution  is  exhibited  in  Fig,  15 
at  *-32.  the  location  of  t lie  density  minimum  in  the  * 
profile.  Again,  a  double  log  plot  is  given  so  that  both 
positive  and  negative  density  excesses  ran  be  conve¬ 
nient  lv  represented.  This  time  a  density  pile  up  is  ob¬ 
served  because  the  fluid  has  been  blown  away  from  the 
axis  so  fast  Miat  a  compression  wave  Is  generated.  The 
zero  in  the  ^.’apli  is  at  0.  36);,.  right  at  the  half-widths 
of  the  thermal  and  intensity  distributions.  Thus  inside 
the  thermal  half-width  the  density  is  depressed;  outside 
the  fluid  has  piled  up.  The  radial  density  profile  shown 
In  Fig.  15  Is  in  excellent  agreement  with  that  predicted 
In  (27).  For  example,  the  zero  observed  at  r/>0* 0  36 
Is  predicted  to  occur  at 

)7rn=t/vT=  0.354. 

oimilarly,  the  location  of  the  peak  observed  at  r/r„ 

*0.  51  is  predicted  to  occur  at 

>A=i. 

Similarly,  the  ratio  of  peak  height  to  valley  depth  is  also 
correctly  predicted.  As  a  matter  of  fact,  when  one  takes 
the  trouble  to  evaluate  the  function  /•(*,/)  appearing  in 
(27),  he  finds  precise  agreement  between  (27)  and  Fig. 
15.  Thus,  iioth  (27 >  and  (28)  are  in  complete  quantitative 
agreement  with  the  result  of  the  computer  calculation. 

Figures  16  and  17  exhibit  the  parade  of  effects,  illus¬ 
trating  graphically  the  various  delays,  pulse  shapes, 
and  widths.  Phvsicaily  the  delays  make  sense.  First  tiie 
beam  blasts  through,  heating  the  fluid  as  it  passes.  As 
explained  above,  I  ho  temperature  maximum  moves  at 
i  and,  thus,  behind  the  laser  peak.  As  tills  tempera¬ 
ture  wave  passes  along,  the  fluid  picks  up  kinetic  energy 
and  the  flow  velocities  increase.  The  center  of  this  ef¬ 


fect  trails  the  heat  wave,  allowing  time  for  the  fluid  to 
res|Kind.  Then,  as  the  fluid  begins  to  flow  away  from  the 
propagating  center  of  velocity,  density  deficits  ire  left 
In  the  wake  and  corresponding  radial  compression 
waves  set  out  from  the  beam  axis. 

VI.  SUMMARY  OF  RESULTS 

In  (xirtlons  of  the  pulse  where  the  intensity  Is  small, 
there  is  very  little  interaction  with  the  fluid  and  these 
portions  move  without  appreciable  distortion.  The  peak 
of  the  pulse,  however,  interacts  fairly  strongly  with  the 
fluid  and  the  peak  Is  delayed  relative  to  the  center  of  the 
pulse.  A  parade  of  effeits  ensues;  tiie  center  and  edges 
of  tiie  pulse  are  followed  by  the  peak,  which,  in  turn,  is 
followed  consecutively  by  t lie  thermal  wave,  the  center 
of  velocity,  and  the  density  waves,  Tiie  front  edge  of  the 
pulse  propagates  without  appreciable  distortion,  but 
streng  phase  oscillations  are  set  up  near  the  peak  and 
rapidly  overtake  the  undistorted  front  section  indicating 
that  soon  the  entire  beam  will  be  distorted  to  sumo  de¬ 
gree.  The  strongest  instabilities  predicted  In  the  lin¬ 
earized  analysis  of  Ref.  1  could  not  appear  in  the  com¬ 
puter  solution  because  they  are  general ed  by  ripple  with 
wavelength  an  order  of  magnitude  smaller  than  the  mesh 
size  used.  The  success  of  the  analytic  analysis  present¬ 
ed  here  suggests,  however,  that  such  instabilities  are 
not  important  for  the  time  interval  considered. 

There  is  very  little  hope  of  obtaining  computer  solutions 
of  the  laser-fluid  equations  except  in  the  tightly  limited 
regime  reported  here,  unless  a  different  ralcnlatlonal 
procedure  can  lie  devised.  Since,  however,  it  appears 
that  a  certain  amount  of  analytical  headway  has  been 
made,  there  is  reason  to  believe  that,  with  appropriate 
combination  of  analytical  and  computer  methods,  the 
tream  can  be  followed  for  considerably  longer  periods  of 
time.  Kffort  Is  currently  being  directed  toward  this 
objective. 
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